Statistik och dataanalys I, 15 hp

Inlämningsuppgift 1

Authors

Johanna Selö

Minhui Zhong

Publicerad

October 3, 2023

Warning

Den här inlämningsuppgiften förutsätter att följande paket finns installerade:

  • mosaic

  • dplyr

  • geosphere

  • leaflet

Paket kan installeras via kommandot install.packages('packagename'), där 'packagename' är namnet på paketet, t.ex 'mosaic'.

Introduktion

I den första inlämningsuppgiften ska ni självständigt i grupper om tre analysera ett dataset i programmeringsspråket R. Till skillnad från datorlaborationerna finns det minimalt med kodexempel. Datorlaborationerna går igenom de flesta momenten som behandlas i inlämningsuppgiften, så se till att göra klart dessa innan.


Instruktioner

I denna inlämningsuppgift ska ni analysera ett datamaterial som beskriver ca 500 distrikt i Boston år 1970. Datasetet är en modifierad version1 av originaldata2 som användes i en studie3 där författarna predikterade medianhuspriser i olika distrikt med hjälp av en rad förklaringsvariabler.

Följande variabler finns i datasetet boston_census_data.Rdata (ladda ner), som innehåller 480 observationer. Notera att en observation motsvarar ett distrikt:

  • town: Stadsdel.
  • longitude: Longitud koordinat.
  • latitude: Latitud koordinat.
  • median_home_value: Medianhuspriset (enhet 1K USD).
  • crime_rate: Brott (per 1000 invånare).
  • zoned_25k_p: Andel av stadsdelens bostadsmark ämnad för marklotter större än 25000 kvadratfot.
  • indust_p: Andel tunnland ägd av företag utanför detaljhandel.
  • borders_charles: Charles River dummy variabel (= 1 om området gränsar till floden, 0 annars).
  • NOx: Koncentration av kväveoxider (andelar per 10 miljon).
  • n_rooms_avg: Genomsnitt antal rum i ägda bostäder.
  • before_1940_p: Andel ägda bostäder byggda före 1940.
  • employ_dist: Viktat avstånd till fem arbetsförmedlingscentra i Boston.
  • radial_access: Index som mäter tillgång till stadsmotorvägar.
  • tax_rate: Fastighetsskatt per 10000 USD.
  • pupil_teacher_ratio: Lärartäthet mätt som elev per lärare.
  • lower_stat_pct: Procentandel med låg socioekonomisk status i termer av utbildning eller arbete.
  • dist_fenway_park: Avstånd till stadion Fenway Park.

Inlämningsuppgiften ska lämnas in i form av ett html-dokument genererat av Quarto. Kontrollera att ni inte får några felmeddelande när du genererar HTML-dokumentet. Läs sedan igenom HTML-dokumentet noggrant innan ni lämnar in det. Använd tydliga figurer och namnge axlarna med tydliga variabelnamn. Glöm inte att skriva era namn i början av dokumentet, där det nu står Namn 1, Namn 2 och Namn 3.

0. Ladda in data

💪 Uppgift 0.1

Ladda in datasetet Boston_census_data med följande kod.

Uppgift 0.1 - Svar
# Write your code here
load(file = url("https://github.com/StatisticsSU/SDA1/blob/main/assignments/assignment1/Boston_census_data.RData?raw=true")) 

1. Kriminalitet i Boston

I detta avsnitt ska ni analysera kriminaliteten i Boston med hjälp av variabeln crime_rate.

💪 Uppgift 1.1

Vad kan man generellt säga om kriminaliteten i Boston 1970? Använd lämpliga figurer och mått för att ge en beskrivning.

Uppgift 1.1 - Svar

Variabeln crime_rate mäter brott per 1000 invånare. Vi börjar med att använda kommandot “favstats” för att få en bra överblick över variabeln. Favstats säger oss att det minsta värdet är 0.00632 och 88.9762 är det högsta, det ger en varians på 88,96988 vilket i sammanhanget är väldigt mycket och säger oss att brottsligheten skiljer sig väldigt mycket åt mellan olika distrikt. Medelvärdet är 3.66526 med en standardavvikelse på 8.746 medan Medianen enbart är 0.253715 med en IQR på 3.599. I och med att datasetet är snedvridet så anser vi att medianen är ett bättre mått för denna variabel.

Genom att använda en boxplot bekräftar det tanken på att datasetet är väldigt snedvridet då vi ser väldigt många outliers. Vi kan även se om vi gör ett histogram istället att majoriteten av distrikten, över 80%, har en crime_rate på mellan 0 och 10 brott per tusen invånare, vilket med andra ord betyder att de flesta distrikten har en relativt låg crime_rate medan några få distrikt har en väldigt hög crime_rate.

# Write your code here
suppressMessages(library(mosaic))
suppressMessages(library(dplyr))
suppressMessages(library(geosphere))
suppressMessages(library(leaflet))
favstats(~crime_rate, data = Boston_census_data)
     min        Q1   median       Q3     max    mean       sd   n missing
 0.00632 0.0829725 0.253715 3.681942 88.9762 3.66526 8.746156 480       0
bwplot(~crime_rate, data=Boston_census_data)

histogram(~crime_rate, data = Boston_census_data, breaks = 10)

💪 Uppgift 1.2

Distrikten tillhör olika stadsdelar som anges i den kategoriska variabeln town? Det finns 88 olika sådana stadsdelar (towns).

Skiljer sig brottsligheten åt mellan de olika stadsdelarna? Undersök stadsdelarna Boston East Boston, Boston Downtown, Cambridge, samt ytterligare en stadsdel som ni själva väljer. Besvara frågan med hjälp av lämpligt valda figurer och statistiska mått.

Tip

Innan ni påbörjar analysen, skapa en ny data frame som enbart innehåller de stadsdelar som ni vill jämföra. Det kan göras till exempel med funktionen filter().

Uppgift 1.2 - Svar

Efter att vi har skapat vårt filter som enbart innehåller de utvalda stadsdelarna, n=63 istället för 480, kan precis som på förra deluppgiften börja med kommandot favstats. Där ser vi att variansen fortfarande är väldigt stor 51.1358 - 0.29819 = 50.83761 samt att medianen och medelvärdet skiljer sig åt, median = 2.3139 och medelvärdet = 6.261704. Det säger oss att fördelningen fortfarande är skev, både medelvärdet och medianen är högre här än om vi tittar på alla observationer.

Tittar vi på histogrammet ser vi tydligt att den är högerskev med en markant outlier med en crime_rate på 50 per 1000 invånare. Fördelningen är unimodal med ett typvärde på en crime_rate mellan 0 och 10 per 1000 invånare.

Avslutningsvis gör vi en korstabell, där vi ser att om man ställer stadsdelarna mot varandra, att majoriteten av brotten sker I Cambridge (41,27%).

# Write your code here
Boston_census_new <- filter(Boston_census_data, town == "Boston East Boston" |town == "Boston Downtown" |town == "Cambridge" |town == "Newton")
favstats(~crime_rate, data = Boston_census_new)
     min      Q1 median      Q3     max     mean       sd  n missing
 0.29819 0.61913 2.3139 8.07211 51.1358 6.261704 9.129887 63       0
histogram(~crime_rate, data = Boston_census_new)

tally(~town + crime_rate, data = Boston_census_new, format = "percent", margins = TRUE)
                    crime_rate
town                    0.29819    0.31533    0.33045    0.33147    0.35809
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            0.000000   0.000000   0.000000   0.000000   0.000000
  Newton               1.587302   1.587302   1.587302   1.587302   1.587302
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    0.38214    0.40771    0.41238    0.44178    0.44791
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            0.000000   0.000000   0.000000   0.000000   0.000000
  Newton               1.587302   1.587302   1.587302   1.587302   1.587302
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    0.46296    0.51183    0.52058      0.537    0.57529
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            0.000000   0.000000   0.000000   0.000000   0.000000
  Newton               1.587302   1.587302   1.587302   1.587302   1.587302
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                     0.6147    0.62356    1.12658    1.20742    1.22358
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            0.000000   0.000000   1.587302   1.587302   1.587302
  Newton               1.587302   1.587302   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    1.27346    1.34284    1.41385    1.42502    1.49632
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            1.587302   1.587302   1.587302   1.587302   1.587302
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                     1.6566    1.80028    2.14918    2.15505    2.24236
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            1.587302   1.587302   1.587302   1.587302   1.587302
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                     2.3004     2.3139    2.33099    2.36862    2.37934
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            1.587302   1.587302   1.587302   1.587302   1.587302
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    2.44668    2.44953    2.73397    2.77974      2.924
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   0.000000   0.000000
  Cambridge            1.587302   1.587302   1.587302   1.587302   1.587302
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    3.32105    3.53501     4.0974    5.29305    6.96215
  Boston Downtown      0.000000   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   0.000000   0.000000   1.587302   1.587302
  Cambridge            1.587302   1.587302   1.587302   0.000000   0.000000
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    7.40389    7.99248    8.15174    9.18702    11.5779
  Boston Downtown      1.587302   0.000000   0.000000   0.000000   0.000000
  Boston East Boston   0.000000   1.587302   1.587302   1.587302   1.587302
  Cambridge            0.000000   0.000000   0.000000   0.000000   0.000000
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    11.9511    14.0507    14.3337    14.4383    15.8744
  Boston Downtown      1.587302   1.587302   0.000000   1.587302   0.000000
  Boston East Boston   0.000000   0.000000   1.587302   0.000000   1.587302
  Cambridge            0.000000   0.000000   0.000000   0.000000   0.000000
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    16.8118     18.811    20.0849    20.7162    22.5971
  Boston Downtown      0.000000   1.587302   0.000000   1.587302   0.000000
  Boston East Boston   1.587302   0.000000   1.587302   0.000000   1.587302
  Cambridge            0.000000   0.000000   0.000000   0.000000   0.000000
  Newton               0.000000   0.000000   0.000000   0.000000   0.000000
  Total                1.587302   1.587302   1.587302   1.587302   1.587302
                    crime_rate
town                    24.3938    28.6558    51.1358      Total
  Boston Downtown      0.000000   1.587302   1.587302  12.698413
  Boston East Boston   1.587302   0.000000   0.000000  19.047619
  Cambridge            0.000000   0.000000   0.000000  41.269841
  Newton               0.000000   0.000000   0.000000  26.984127
  Total                1.587302   1.587302   1.587302 100.000000

💪 Uppgift 1.3

Vilka två variabler i datasetet Boston_census_data korrelerar mest med brottslighet? Beskriv sambandet mellan brottslighet och var och en av dessa två variabler.

Tip

Kom ihåg att korrelation mäter det linjära sambandet mellan numeriska variabler.

Uppgift 1.3 - Svar

Vi skapar ett nytt dataset som enbart innehåller numeriska variabler som vi döper till “Boston_census_data_15_variables”. Sen för att hitta vilka två variabler som korrelerar mest med brottslighet gör vi en correlation matrix, denna är dock svåröversiktlig så då är det bättre att göra en corrplot, där man tydligare ser vad som korrelerar med vad. Detta kan man senare dubbelkolla med correlation matrixen, ifall vissa pluppar är snarlika varandra.

Vi hittar 3 variabler som ser ut att korrelera mest med crime_rate: radial_access, median_home_value och tax_rate. Vi använder funktionen cor() för att se vilka två som har högst samband. radial_access: 0.62, median_home_value: -0.45 och tax_rate: 0.58.

Radial_access och crime_rate har ett relativt starkt positivt samband på 0.62, det innebär ju större tillgång till stadsmotorvägar desto högre crime_rate.

tax_rate och crime_rate har också ett postivt samband på 0.58, det innebär att ju högre tax_rate desto mer brott begås per 1000 invånare.

# Write your code here
Boston_census_data_15_variables <- Boston_census_data[, c("longitude","latitude","median_home_value","crime_rate","zoned_25k_p","indust_p","NOx","n_rooms_avg","before_1940_p","employ_dist","radial_access","tax_rate","pupil_teacher_ratio","lower_stat_pct","dist_fenway_park")]
correlation_matrix_Boston <- cor(Boston_census_data_15_variables)
round(correlation_matrix_Boston, 3)
                    longitude latitude median_home_value crime_rate zoned_25k_p
longitude               1.000    0.140            -0.326      0.056      -0.168
latitude                0.140    1.000             0.015     -0.092      -0.123
median_home_value      -0.326    0.015             1.000     -0.450       0.395
crime_rate              0.056   -0.092            -0.450      1.000      -0.196
zoned_25k_p            -0.168   -0.123             0.395     -0.196       1.000
indust_p                0.035   -0.063            -0.600      0.403      -0.522
NOx                     0.138   -0.085            -0.524      0.415      -0.506
n_rooms_avg            -0.227   -0.069             0.680     -0.213       0.299
before_1940_p           0.183    0.060            -0.488      0.350      -0.559
employ_dist             0.011   -0.069             0.368     -0.378       0.672
radial_access           0.015   -0.226            -0.479      0.623      -0.303
tax_rate                0.031   -0.189            -0.579      0.580      -0.300
pupil_teacher_ratio     0.295   -0.002            -0.515      0.283      -0.375
lower_stat_pct          0.172    0.034            -0.765      0.464      -0.424
dist_fenway_park        0.435   -0.200            -0.003     -0.123       0.400
                    indust_p    NOx n_rooms_avg before_1940_p employ_dist
longitude              0.035  0.138      -0.227         0.183       0.011
latitude              -0.063 -0.085      -0.069         0.060      -0.069
median_home_value     -0.600 -0.524       0.680        -0.488       0.368
crime_rate             0.403  0.415      -0.213         0.350      -0.378
zoned_25k_p           -0.522 -0.506       0.299        -0.559       0.672
indust_p               1.000  0.761      -0.407         0.632      -0.706
NOx                    0.761  1.000      -0.315         0.725      -0.766
n_rooms_avg           -0.407 -0.315       1.000        -0.256       0.242
before_1940_p          0.632  0.725      -0.256         1.000      -0.744
employ_dist           -0.706 -0.766       0.242        -0.744       1.000
radial_access          0.592  0.608      -0.190         0.447      -0.486
tax_rate               0.717  0.666      -0.278         0.499      -0.529
pupil_teacher_ratio    0.382  0.178      -0.289         0.264      -0.237
lower_stat_pct         0.640  0.615      -0.607         0.633      -0.544
dist_fenway_park      -0.394 -0.325       0.020        -0.304       0.686
                    radial_access tax_rate pupil_teacher_ratio lower_stat_pct
longitude                   0.015    0.031               0.295          0.172
latitude                   -0.226   -0.189              -0.002          0.034
median_home_value          -0.479   -0.579              -0.515         -0.765
crime_rate                  0.623    0.580               0.283          0.464
zoned_25k_p                -0.303   -0.300              -0.375         -0.424
indust_p                    0.592    0.717               0.382          0.640
NOx                         0.608    0.666               0.178          0.615
n_rooms_avg                -0.190   -0.278              -0.289         -0.607
before_1940_p               0.447    0.499               0.264          0.633
employ_dist                -0.486   -0.529              -0.237         -0.544
radial_access               1.000    0.910               0.453          0.513
tax_rate                    0.910    1.000               0.455          0.572
pupil_teacher_ratio         0.453    0.455               1.000          0.364
lower_stat_pct              0.513    0.572               0.364          1.000
dist_fenway_park           -0.206   -0.233               0.043         -0.161
                    dist_fenway_park
longitude                      0.435
latitude                      -0.200
median_home_value             -0.003
crime_rate                    -0.123
zoned_25k_p                    0.400
indust_p                      -0.394
NOx                           -0.325
n_rooms_avg                    0.020
before_1940_p                 -0.304
employ_dist                    0.686
radial_access                 -0.206
tax_rate                      -0.233
pupil_teacher_ratio            0.043
lower_stat_pct                -0.161
dist_fenway_park               1.000
suppressMessages(library(corrplot))
corrplot(correlation_matrix_Boston)

cor(radial_access ~ crime_rate, data = Boston_census_data)
[1] 0.623462
cor(median_home_value ~ crime_rate, data = Boston_census_data)
[1] -0.4496238
cor(tax_rate ~ crime_rate, data = Boston_census_data)
[1] 0.5802444

2. Fastighetsskatt i Boston

I detta avsnitt ska ni analysera fastighetsskatten i Boston med hjälp av variabeln tax_rate.

💪 Uppgift 2.1

Vad kan man generellt säga om fastighetsskatten i distrikten? Använd lämpliga figurer och mått för att beskriva fördelningen.

Uppgift 2.1 - Svar

Fastighetsskatten varierar mellan 187 och 711 per 1000 USD, vilket ger en varians på 524. Medianen och Medelvärdet skiljer sig inte lika mycket som det gjorde på brottsligheten, utan ligger här på 330 respektive 409.33. Standardavvikelsen är 168.5777 och IQR är 385.

Om vi gör ett histogram ser vi att det är en bimodal fördelning samt att det ser ut att vara två olika grupper i och med att det är ett glapp mellan staplarna. Typvärderna ligger runt 300 och 650. Antagligen finns det ett samband med en annan variabel som avgör om man betalar lägre eller högre fastighetsskatt.

# Write your code here
favstats(~tax_rate, data=Boston_census_data)
 min     Q1 median  Q3 max     mean       sd   n missing
 187 280.75    330 666 711 409.3271 168.5777 480       0
histogram(~ tax_rate, data = Boston_census_data, col = "navyblue", type="count")

💪 Uppgift 2.2

Låt oss skapa en ny variabel cat_tax som anger om ett distrikt betalar låg (low), medel (medium), eller hög (high) fastighetsskatt. Vi definerar skattekategorierna enligt

  • low: tax_rate \(\leq\) 250,
  • medium: 250 \(<\) tax_rate \(\leq\) 400,
  • high: tax_rate \(>\) 400.

Följande kod skapar och lägger till variabeln cat_tax i Boston_census_data

Boston_census_data$cat_tax <- cut(Boston_census_data$tax_rate, 
              breaks=c(0, 250, 400, 800),
              labels=c('Low', 'Medium', 'High'))

Finns det ett samband mellan vilken skattekategori ett distrikt tillhör och dess angränsning till Charles River? Förklara med hjälp av lämplig tabell och figur.

Uppgift 2.2 - Svar

Som man ser i stapeldiagrammet skiljer sig skattenivåerna inte så mycket åt mellan om man gränsar till Charles River eller inte. Det man kan se är att det är mer troligt att man har men medium-tax rate om man gränsar till Charles River än om man inte gör det.

Fördelningen mellan skattenivåer om man inte gränsar till Charles River eller om man gör det kan vi se i korstabellen. Det är större sannolikhet att man inte gränsar till Charles River om man har låg (13,3% mot 10,3%) eller hög (39,69% mot 34,48%) skattekategori.

# Write your code here
tally(cat_tax ~ borders_charles, data = Boston_census_data, margins=TRUE, format="percent")
        borders_charles
cat_tax          0         1
  Low     13.30377  10.34483
  Medium  47.00665  55.17241
  High    39.68958  34.48276
  Total  100.00000 100.00000
bargraph(~cat_tax, groups=borders_charles, data=Boston_census_data, type="percent", main = "Skattekategorier baserat på angräsning till Charles River")

💪 Uppgift 2.3

Hur många procent av alla distrikt i vår data ligger i angränsning till Charles River och tillhör en hög skattekategori? Hur stor andel av distrikten med hög skatt ligger inte i angränsning till Charles River?

Uppgift 2.3 - Svar

Nedan ser vi först en marginalfördelningstabell. Där kan vi se att 2,08% av alla distrikt i vår data angränsar till River Charles (1) och tillhör en hög skattekategori (High).

Vi ser sedan en betingad fördelningstabell, som är uppdelad på skattekategorierna. I den kan vi se att av de distrikt som har en hög skattekategori är det 94,71% som inte ligger i angränsning till River Charles (0).

# Write your code here
tally(~ cat_tax + borders_charles, data=Boston_census_data, format="percent", margins=T)
        borders_charles
cat_tax           0          1      Total
  Low     12.500000   0.625000  13.125000
  Medium  44.166667   3.333333  47.500000
  High    37.291667   2.083333  39.375000
  Total   93.958333   6.041667 100.000000
tally(borders_charles ~ cat_tax, data = Boston_census_data, margins=TRUE, format="percent")
               cat_tax
borders_charles        Low     Medium       High
          0      95.238095  92.982456  94.708995
          1       4.761905   7.017544   5.291005
          Total 100.000000 100.000000 100.000000

💪 Uppgift 2.4

Vilka två variabler i datasetet Boston_census_data korrelerar starkast med tax_rate? Beskriv det parvisa sambandet mellan tax_rate och var och en av dessa två variabler. Vad kan vi säga om kausalitet för vart och ett av sambanden?

Tip

Kom ihåg att korrelation är ett mått på linjära samband mellan numeriska variabler.

Uppgift 2.4 - Svar

Skriv svaret här.

Vi skapar ett nytt dataset som enbart innehåller numeriska variabler som vi döper till “Boston_census_data_15_variables”. Sen för att hitta vilka två variabler som korrelerar mest med fastighetsskatt gör vi en correlation matrix, denna är dock svåröversiktlig så då är det bättre att göra en corrplot, där man tydligare ser vad som korrelerar med vad. Detta kan man senare dubbelkolla med correlation matrixen, ifall vissa pluppar är snarlika varandra.

Man ser ganska tydligt att det är två variabler som verkar korrelera mest med tax_rate och det är radial_access och indust_p. Vi använder cor()-funktionen för att se korrelationerna.

Radial_access och tax_rate har ett väldigt starkt positivt samband på 0.91. Det verkar finnas ett tydligt samband mellan större tillgång till stadsmotorvägar och distriktens fastighetsskatt.

Även indust_p och tax_rate har ett starkt positivt samband på 0.72. indust_p står för andel tunnland ägd av företag utanför detaljhandeln. Det kan man tolka att ju mer tunnland som ägs av företag utanför detaljhandeln, desto högre tax_rate.

Vi kan inte säga att någon av dessa är kausala samband, vi vet inte med säkerhet att det är t ex tillgången till stadsmotorvägar som påverkar fastighetsskatten, utan det kan alltid finnas en dold variabel utanför vårt dataset som påverkar både, och bara får det att se ut som att det finns ett samband.

# Write your code here
Boston_census_data_15_variables <- Boston_census_data[, c("longitude","latitude","median_home_value","crime_rate","zoned_25k_p","indust_p","NOx","n_rooms_avg","before_1940_p","employ_dist","radial_access","tax_rate","pupil_teacher_ratio","lower_stat_pct","dist_fenway_park")]
suppressMessages(library(corrplot))
corrplot(correlation_matrix_Boston)

cor(radial_access ~ tax_rate, data = Boston_census_data)
[1] 0.9099416
cor(indust_p ~ tax_rate, data = Boston_census_data)
[1] 0.7171041

3. Avstånd till Fenway park

I detta avsnitt ska ni undersöka variabeln dist_fenway_park, som mäter avståndet mellan ett distrikt och Fenway park (stadion där basebollslaget Boston Red Sox spelar sina hemmamatcher).

Vi kan visualisera Fenway park och distrikten på en karta med hjälp av R-paketet leaflet. Följande kod visar platsen för Fenway park och distrikten för observationerna 30 och 45.

library(leaflet) # Install if not available
fenway_park_lat_long <- c(42.346462, -71.097250) # latitude and longitude for Fenway_park
Boston_map <- leaflet() %>% 
  addTiles() %>%
  addMarkers(lat = fenway_park_lat_long[1], lng = fenway_park_lat_long[2], popup="Fenway park") %>%
  addMarkers(lat = Boston_census_data$latitude[30], lng = Boston_census_data$longitude[30], popup="Observation 30") %>%
  addMarkers(lat = Boston_census_data$latitude[45], lng = Boston_census_data$longitude[45], popup="Observation 45") 

Boston_map # Show interactive map

💪 Uppgift 3.1

Vilket distrikt i vår data har längst respektive kortast avstånd till Fenway park? Markera ut dessa distrikt i en interaktiv karta tillsammans med Fenway park.

Uppgift 3.1 - Svar

Det första vi måste göra är att få fram vilka distrikt som har längst respektive kortast avstånd till Fenway park. Det får vi fram genom att sortera datan på “dist_fenway_park”-variablen. Om vi använder funktionen head() kan vi se distrikten på först raden. Wilmington har kortast avstånd på 887.9 medan Marshfield har längst på 33638.4.

För att få ut dessa på kartan, betyder vi ut till våra nya variabler för den sorterade datan och väljer den första raden i respektive.

dist_fenway_sorted <- Boston_census_data %>% arrange(dist_fenway_park)
dist_fenway_sorted_hightolow <- Boston_census_data %>% arrange(desc(dist_fenway_park))

library(leaflet) 
fenway_park_lat_long <- c(42.346462, -71.097250) # latitude and longitude for Fenway_park
Boston_map <- leaflet() %>% 
  addTiles() %>%
  addMarkers(lat = fenway_park_lat_long[1], lng = fenway_park_lat_long[2], popup="Fenway park") %>%
  addMarkers(lat = dist_fenway_sorted$latitude[1], lng = dist_fenway_sorted$longitude[1], popup="kortast") %>%
  addMarkers(lat = dist_fenway_sorted_hightolow$latitude[1], lng = dist_fenway_sorted_hightolow$longitude[1], popup="Observation längst") 
Boston_map

💪 Uppgift 3.2

Finns det ett samband mellan dist_fenway_park och crime_rate?

Uppgift 3.2 - Svar

Det finns ett svagt negativt samband mellan avståndet till Fenway park och crime_rate. Korrelationen antyder att närmare Fenway park (lägre avstånd) desto mer brottslighet (högre crime_rate). Dock är sambandet svagt så det kan likaväl bero på slumpen.

# Write your code here
cor(dist_fenway_park ~ crime_rate, data = Boston_census_data)
[1] -0.1228651

4. Enkel linjär regression

I detta avsnitt ska ni anpassa och tolka enkla linjära regressionsmodeller.

💪 Uppgift 4.1

Anpassa en linjär regression med responsvariabeln NOx och den förklarande variabeln employ_dist. Rita den anpassade regressionslinjen tillsammans med data i en lämplig figur. Beskriv resultaten och tolka modellen. Utför en modellvalidering via en residualanalys och kommentera modellens lämplighet. Om modellen inte anses lämplig, vilka antaganden har inte varit uppfyllda?

Uppgift 4.1 - Svar

Skriv svaret här.

# Write your code here
lm_NOx_vs_employdist <- lm(NOx ~ employ_dist, data = Boston_census_data)

summary(lm_NOx_vs_employdist)

Call:
lm(formula = NOx ~ employ_dist, data = Boston_census_data)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.12358 -0.05352 -0.01220  0.04348  0.22868 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.717614   0.007109  100.95   <2e-16 ***
employ_dist -0.042641   0.001639  -26.02   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.07526 on 478 degrees of freedom
Multiple R-squared:  0.5861,    Adjusted R-squared:  0.5852 
F-statistic: 676.9 on 1 and 478 DF,  p-value: < 2.2e-16
plot(NOx ~ employ_dist, data = Boston_census_data, col = "cornflowerblue", ylim = c(0, 1))
y_hat <- predict(lm_NOx_vs_employdist)
head(y_hat)
        1         2         3         4         5         6 
0.6337437 0.5124177 0.6522413 0.6235483 0.6111782 0.4921591 
lines(Boston_census_data$employ_dist, y_hat, type = "p", col = "lightcoral")
abline(lm_NOx_vs_employdist, col = "lightcoral")
legend(x = "topleft", pch = c(1, 1, NA), lty = c(NA, NA, 1), col = c("cornflowerblue", "lightcoral", "lightcoral"), legend=c("Data", "Predicted", "Fitted line"))

# residualer
resid <- residuals(lm_NOx_vs_employdist)
head(resid)
           1            2            3            4            5            6 
-0.009743748  0.002582251  0.047758653 -0.039548322  0.101821788 -0.079159083 
# residualanalys
Boston_census_data$res <- resid(lm_NOx_vs_employdist)
Boston_census_data$y_hatt <- fitted(lm_NOx_vs_employdist)
plot(Boston_census_data$res ~ Boston_census_data$y_hatt, ylab="Resid", xlab="y-hatt", main="Residplot")
abline(h=0)

qqnorm(resid, col = "cornflowerblue") # Create normal probability plot for residuals
qqline(resid, col = "red") # Add a straight line to normal probability plot 

# Den är inte så slumpmässig, inte normalfördelad.

💪 Uppgift 4.2

Använd modellen i Uppgift 4.1 för att prediktera koncentration av kväveoxider för observation 10, där employ_dist=10.5857. Beräkna vad residualen blir för denna observation.

Uppgift 4.2 - Svar

Skriv svaret här.

# Write your code here
new_x <- data.frame(employ_dist = c(10.5857))
predict(lm_NOx_vs_employdist, newdata = new_x)
        1 
0.2662308 
# Beräkna residualen

💪 Uppgift 4.3

Transformera variablerna i Uppgift 4.1 (avgör själv vilken eller vilka av de två som behöver transformeras). Ett förslag är att använda Tukeys cirkel för att hitta lämpliga transformationer. Anpassa en ny linjär regression med de transformerade variablerna. Utför en modellvalidering (efter transformation) via en residualanalys och kommentera modellens lämplighet jämfört med modellen i Uppgift 4.1.

Uppgift 4.3 - Svar

Skriv svaret här.

# Write your code here

plot(NOx ~ employ_dist, data = Boston_census_data,col = "cornflowerblue") #otransformerad

# Testar olika transformationer
plot(sqrt(NOx) ~ sqrt(employ_dist), data = Boston_census_data, col = "cornflowerblue")

plot(sqrt(NOx) ~ employ_dist, data = Boston_census_data, col = "cornflowerblue")

plot(NOx ~ sqrt(employ_dist), data = Boston_census_data, col = "cornflowerblue")

plot(log(NOx) ~ log(employ_dist), data = Boston_census_data, col = "pink")

plot(NOx ~ log(employ_dist), data = Boston_census_data, col = "pink")

plot(log(NOx) ~ employ_dist, data = Boston_census_data, col = "pink")

plot(sqrt(NOx) ~ log(employ_dist), data = Boston_census_data, col = "coral")

plot(log(NOx) ~ sqrt(employ_dist), data = Boston_census_data, col = "coral")

# Vi väljer att transformera både x och y till log(x) och log(y) för det ger den rätaste linjen.

lm_logNOx_vs_logemploy_dist <- lm(log(NOx) ~ log(employ_dist), data = Boston_census_data)
summary(lm_logNOx_vs_logemploy_dist)

Call:
lm(formula = log(NOx) ~ log(employ_dist), data = Boston_census_data)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.19152 -0.07828 -0.01370  0.06108  0.28496 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)      -0.218492   0.011524  -18.96   <2e-16 ***
log(employ_dist) -0.327303   0.008832  -37.06   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.103 on 478 degrees of freedom
Multiple R-squared:  0.7418,    Adjusted R-squared:  0.7413 
F-statistic:  1373 on 1 and 478 DF,  p-value: < 2.2e-16
plot(log(NOx) ~ log(employ_dist), data = Boston_census_data, col = "cornflowerblue")
log_y_hat <- predict(lm_logNOx_vs_logemploy_dist)
head(log_y_hat)
         1          2          3          4          5          6 
-0.4398989 -0.7327355 -0.3583458 -0.4774478 -0.5178857 -0.7635522 
plot(log(NOx) ~ log(employ_dist), data = Boston_census_data, col = "cornflowerblue")
lines(Boston_census_data$employ_dist, log_y_hat, type ="p", col = "lightcoral")
abline(lm_logNOx_vs_logemploy_dist, col = "lightcoral")
legend(x = "topleft", pch = c(0,5, 0,5, NA), lty = c(NA, NA, 0,5), col = c("cornflowerblue", "lightcoral", "lightcoral"), legend=c("Data", "Predicted", "Fitted line"))

# Residualanalys
# residualer
resid <- residuals(lm_logNOx_vs_logemploy_dist)
head(resid)
           1            2            3            4            5            6 
-0.031706019  0.069147161  0.001670863 -0.060406468  0.179611805 -0.120755463 
# residualanalys
Boston_census_data$res <- resid(lm_logNOx_vs_logemploy_dist)
Boston_census_data$y_hatt <- fitted(lm_logNOx_vs_logemploy_dist)
plot(Boston_census_data$res ~ Boston_census_data$y_hatt, ylab="Resid", xlab="y-hatt", main="Residplot")
abline(h=0)

qqnorm(resid, col = "cornflowerblue") # Create normal probability plot for residuals
qqline(resid, col = "red") # Add a straight line to normal probability plot 

💪 Uppgift 4.4

Plotta den anpassade regressionen från 4.3 i icke-transformerad skala tillsammans med observationerna (också i icke-transformerad skala) i en lämplig figur.

Uppgift 4.4 - Svar

Skriv svaret här.

# Write your code here

logy_hat <- predict(lm_logNOx_vs_logemploy_dist) # log scale prediction
y_hat <- exp(logy_hat) # original scale prediction
head(logy_hat)
         1          2          3          4          5          6 
-0.4398989 -0.7327355 -0.3583458 -0.4774478 -0.5178857 -0.7635522 
head(y_hat)
        1         2         3         4         5         6 
0.6441015 0.4805925 0.6988314 0.6203647 0.5957789 0.4660081 
plot(NOx ~ employ_dist, data = Boston_census_data, col = "cornflowerblue") # Data on original scale
lines(Boston_census_data$employ_dist, y_hat, type = "p", col = "lightcoral")
legend(x = "topleft", pch = c(1, 1), col = c("cornflowerblue", "lightcoral"), legend=c("Data", "Predicted"))

💪 Uppgift 4.5

Använd modellen i Uppgift 4.3 för att prediktera koncentration av kväveoxider i icke-transformerad skala för observation 10, där employ_dist=10.5857. Beräkna vad residualen blir för denna observation. Kommentera resultaten jämfört med Uppgift 4.2.

Tip

Tänk på att ta hänsyn till eventuella transformationer!

Uppgift 4.5 - Svar

Skriv svaret här.

# Write your code here
new_x <- data.frame(employ_dist = c(10.5857))
predict(lm_logNOx_vs_logemploy_dist, newdata = new_x)
         1 
-0.9907648 
# Beräkna residualen
resid <- residuals(lm_logNOx_vs_logemploy_dist)

5. Multipel linjär regression

I detta avsnitt ska ni studera multipel linjära regression.

💪 Uppgift 5.1

Anpassa en linjär regression med responsvariabel logaritmerad median_home_value samt förklarande variabler lower_stat_pct och dummy-variabeln borders_charles. Tolka koefficienten för borders_charles.

Uppgift 5.1 - Svar

Skriv svaret här.

Om borders_charles är 1, ökar vår prediktion av y (median_home_value) med 2.37.

# Write your code here

lm_median_home_value_vs_lower_stat_pct_borders_charles <- lm(median_home_value ~ lower_stat_pct + borders_charles, data = Boston_census_data)
summary(lm_median_home_value_vs_lower_stat_pct_borders_charles)

Call:
lm(formula = median_home_value ~ lower_stat_pct + borders_charles, 
    data = Boston_census_data)

Residuals:
    Min      1Q  Median      3Q     Max 
-9.2067 -3.1683 -0.9976  1.6929 21.4235 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)     32.36661    0.48051  67.358   <2e-16 ***
lower_stat_pct  -0.84435    0.03231 -26.130   <2e-16 ***
borders_charles  2.37058    0.95755   2.476   0.0136 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 4.998 on 477 degrees of freedom
Multiple R-squared:  0.5912,    Adjusted R-squared:  0.5895 
F-statistic:   345 on 2 and 477 DF,  p-value: < 2.2e-16
# Göra en regressionslinje

resid_median_multiplereg <- residuals(lm_median_home_value_vs_lower_stat_pct_borders_charles)
head(resid_median_multiplereg)
         1          2          3          4          5          6 
-0.6971156 -5.8226098  0.6705960  1.2821955 -2.1688443 -2.3244440 
plot(lm_median_home_value_vs_lower_stat_pct_borders_charles$fitted.values, resid_median_multiplereg, xlab= "y_hat", ylab='Residuals', col = "cornflowerblue") 

qqnorm(resid_median_multiplereg, col = "cornflowerblue") # Create normal probability plot for residuals
qqline(resid_median_multiplereg, col = "red") # Add a straight line to normal

lmod_boston <- lm(median_home_value ~ lower_stat_pct + borders_charles, data=Boston_census_data)
lmod_boston

Call:
lm(formula = median_home_value ~ lower_stat_pct + borders_charles, 
    data = Boston_census_data)

Coefficients:
    (Intercept)   lower_stat_pct  borders_charles  
        32.3666          -0.8443           2.3706  

💪 Uppgift 5.2

Ni ska nu utforma en modell som predikterar medianhuspriset median_home_value. Ni får endast använda följande förklaringsvariabler:

  • before_1940_p
  • crime_rate
  • radial_access
  • NOx
  • dist_fenway_park

Ni får själva välja hur många av variablerna som ska ingå i modellen. Ni får göra vilka transformationer ni vill av variablerna, inklusive responsvariabeln.

Pröva er fram metodiskt när ni väljer vilka variabler ni inkluderar i modellen, och när ni bestämmer vilka eventuella transformationer ni använder.

När ni utvärderar olika modeller kan ni förslagsvis börja med att jämföra adjusted R-squared. När ni med hjälp av adjusted R-squared har identifierat två eller tre modeller som ser lovande ut kan ni utvärdera dessa modeller ytterligare i ett andra steg.

I det andra steget ska ni utvärdera hur väl modellerna predikterar data som inte använts för att anpassa modellen. Ni kan välja en av två alternativa metoder:

  1. Dela in ert dataset i träningsdata och testdata. Anpassa modellen med hjälp av träningsdata, och utvärdera sedan på testdata. Ni kan exempelvis använda de första 350 observationerna som träningsdata och de sista 130 observationerna som testdata.
  2. Använd korsvalidering. Det är en något mer krävande metod, men också något bättre. Ni kan exempelvis göra korsvalidering med 4 folds (4-fold cross validation). Dela då upp ert dataset i fyra delar (del 1: observationer 1-120, del 2: observationer 121-240, del 3: observationer 241-360, del 4: observationer 361-480).

Sortera inte observationerna i Boston_census_data slumpmässigt. Ordningen är redan slumpmässig.

Tip

Tänk på att ta hänsyn till eventuell transformation av responsvariabeln. Om ni exempelvis har valt transformationen \(\log(y)\) är modellens prediktion av responsvariabeln \(\widehat{\log(y)}\). Ni måste då transformera den till \(\hat y\) i responsvariabelns originalskala med formeln \(\hat{y}=\exp\left(\widehat{\log(y)}\right)\). Sedan kan ni räkna ut residualen \(y - \hat y\).

Uppgift 5.2 - Svar

Skriv svaret här.

# Write your code here

#Model 1: median_home_value vs before_1940_p, model 2: median_home_value vs before_1940p + NOx
model1 <- lm(sqrt(median_home_value) ~ before_1940_p, data=Boston_census_data) #Träna modellen
summary(model1)

Call:
lm(formula = sqrt(median_home_value) ~ before_1940_p, data = Boston_census_data)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.91986 -0.44131 -0.06919  0.34746  2.76484 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)    5.61658    0.08675   64.75   <2e-16 ***
before_1940_p -0.01525    0.00117  -13.04   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.7184 on 478 degrees of freedom
Multiple R-squared:  0.2624,    Adjusted R-squared:  0.2609 
F-statistic: 170.1 on 1 and 478 DF,  p-value: < 2.2e-16
model2 <- lm(sqrt(median_home_value) ~ before_1940_p + NOx, data=Boston_census_data) #Träna modellen
summary(model2) #Se sammanfattande resultat

Call:
lm(formula = sqrt(median_home_value) ~ before_1940_p + NOx, data = Boston_census_data)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.73795 -0.41932 -0.08873  0.34342  2.82617 

Coefficients:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)    6.573798   0.158582  41.453  < 2e-16 ***
before_1940_p -0.006952   0.001618  -4.296 2.11e-05 ***
NOx           -2.748253   0.388636  -7.072 5.48e-12 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.6842 on 477 degrees of freedom
Multiple R-squared:  0.3324,    Adjusted R-squared:  0.3296 
F-statistic: 118.8 on 2 and 477 DF,  p-value: < 2.2e-16
#1 Dela dataset i 350 observationer som träningsdata och 130 observationer som testdata
library(dplyr)
data_randomorder <- Boston_census_data %>% slice_sample(prop=1)
datatrain <- data_randomorder %>% slice(1:350)
datatest <- data_randomorder %>% slice(351:480)

#Utvärdera på testdata
rsquared1 <- summary(model1)$r.squared
y_hatt_test <- predict(model1, newdata=Boston_census_data)
y_test <- Boston_census_data$y
SSE <- sum((y_test - y_hatt_test)^2)
n_test <- 130
MSE <- SSE / n_test
RMSE <- sqrt(MSE)
sprintf("R-squared=%.7f, SSE=%.7f", rsquared1, SSE)
[1] "R-squared=0.2624368, SSE=13025.6687976"
rsquared1 <- summary(model2)$r.squared
y_hatt_test <- predict(model2, newdata=Boston_census_data)
y_test <- Boston_census_data$y
SSE <- sum((y_test - y_hatt_test)^2)
n_test <- 130
MSE <- SSE / n_test
RMSE <- sqrt(MSE)
sprintf("R-squared=%.7f, SSE=%.7f", rsquared1, SSE)
[1] "R-squared=0.3324226, SSE=13063.6016674"
#2 Korsvalidering för Modellen 
n <- 480 # Number of observations
# Fold 1:
obs_index <- c(1:n) # Keeps track of the indices of  the dataset (1, 2, 3, ...., n = 480)
test_fold_index <- obs_index[c(1:120)] # Subsets indices 1:120 (test data fold 1) 
training_fold_index <- obs_index[-c(1:120)] # Takes out the complement
lm_modell1_fold1 <- lm(sqrt(median_home_value) ~ before_1940_p, subset = training_fold_index, data = Boston_census_data) # Estimate fold 1
test_data <- Boston_census_data[test_fold_index, ] # Create test data for fold
sqrty_hat_fold1 <- predict(lm_modell1_fold1, newdata = test_data) # Predict test data in sqrt scale
y_hat_fold1 <- (sqrty_hat_fold1)^2 # Predict test data ordinary scale
SSE_fold1 <- sum((test_data$median_home_value - y_hat_fold1)^2) 
# Fold 2:
test_fold_index <- obs_index[c(121:240)] # Subsets indices 121:240 (test data fold 2) 
training_fold_index <- obs_index[-c(121:240)] # Takes out the complement
lm_modell1_fold2 <- lm(sqrt(median_home_value) ~ before_1940_p, subset = training_fold_index, data = Boston_census_data) # Estimate fold 2
test_data <- Boston_census_data[test_fold_index, ] # Create test data for fold
sqrty_hat_fold2 <- predict(lm_modell1_fold2, newdata = test_data) # Predict test data in sqrt scale
y_hat_fold2 <- (sqrty_hat_fold2)^2 # Predict test data ordinary scale
SSE_fold2 <- sum((test_data$median_home_value - y_hat_fold2)^2)
# Fold 3:
test_fold_index <- obs_index[c(241:360)] # Subsets indices 241:360 (test data fold 3) 
training_fold_index <- obs_index[-c(241:360)] # Takes out the complement
lm_modell1_fold3 <- lm(sqrt(median_home_value) ~ before_1940_p, subset = training_fold_index, data = Boston_census_data) # Estimate fold 3
test_data <- Boston_census_data[test_fold_index, ] # Create test data for fold
sqrty_hat_fold3 <- predict(lm_modell1_fold3, newdata = test_data) # Predict test data in sqrt scale
y_hat_fold3 <- (sqrty_hat_fold3)^2 # Predict test data ordinary scale
SSE_fold3 <- sum((test_data$median_home_value - y_hat_fold3)^2)
# Fold 4:
test_fold_index <- obs_index[c(361:480)] # Subsets indices 361:480 (test data fold 4) 
training_fold_index <- obs_index[-c(361:480)] # Takes out the complement
lm_modell1_fold4 <- lm(sqrt(median_home_value) ~ before_1940_p, subset = training_fold_index, data = Boston_census_data) # Estimate fold 4
test_data <- Boston_census_data[test_fold_index, ] # Create test data for fold
sqrty_hat_fold4 <- predict(lm_modell1_fold4, newdata = test_data) # Predict test data in sqrt scale
y_hat_fold4 <- (sqrty_hat_fold4)^2 # Predict test data ordinary scale
SSE_fold4 <- sum((test_data$median_home_value - y_hat_fold4)^2)

n <- nrow(Boston_census_data)
RMSE = sqrt((SSE_fold1 + SSE_fold2 +SSE_fold3 +SSE_fold4)/n)
RMSE
[1] 6.891041

💪 Uppgift 5.3

Gör en residualanalys av den valda modellen i Uppgift 5.2.

Uppgift 5.3 - Svar

Skriv svaret här.

# Write your code here
# Residualanalys
# residualer
resid <- residuals(model2)
head(resid)
          1           2           3           4           5           6 
-0.03750447 -0.42738065 -0.44764152 -0.64473427 -0.16533605 -0.65877644 
Boston_census_data$res <- resid(model2)
Boston_census_data$y_hatt <- fitted(model2)
plot(Boston_census_data$res ~ Boston_census_data$y_hatt, ylab="Resid", xlab="y-hatt", main="Residplot")
abline(h=0)

qqnorm(resid, col = "cornflowerblue") # Create normal probability plot for residuals
qqline(resid, col = "red") # Add a straight line to normal probability plot

💪 Uppgift 5.4

Använd modellen i Uppgift 5.2 för att prediktera medianhuspriset för observationerna i datasetet Boston_districts_to_predict (ladda ner).

Uppgift 5.4 - Ladda hem data för prediktion

Ladda in dataseten Boston_districts_to_predict med följande kod.

# Write your code here
load(file = url("https://github.com/StatisticsSU/SDA1/blob/main/assignments/assignment1/Boston_districts_to_predict.RData?raw=true")) 

Det här datasetet har endast de förklarande variablerna, dvs inte responsvariabeln. När vi rättar era inlämningsuppgifter kommer vi att jämföra era prediktioner med de faktiska medianpriserna (som vi har tillgång till).

Skriv ut dina prediktioner så att vi enkelt kan se dem när vi rättar.

Tip

Tänk på att ta hänsyn till eventuella transformationer av förklaringsvariablerna!

Uppgift 5.4 - Svar

Skriv svaret här.

# Write your code here

new_x <- data.frame(Boston_districts_to_predict)
predict(model2, newdata = new_x)
       1        2        3        4        5        6        7        8 
5.102926 3.484826 5.341043 4.585608 5.011061 4.214070 4.842807 4.207502 
       9       10 
4.015959 4.800209 
# Transformera tillbaka

Footnotes

  1. Totalundersökningen trunkerade medianhusvärdet till 50K för de censusdistrikten som låg över. Vi har tagit bort dessa censusdistrikt. Vi har också tagit bort variabler som är irrelevanta.↩︎

  2. Harrison Jr, D., & Rubinfeld, D. L. (1978). Hedonic housing prices and the demand for clean air. Journal of Environmental Economics and Management, 5(1), 81-102.↩︎

  3. Pace, R. K., & Gilley, O. W. (1997). Using the spatial configuration of the data to improve estimation. The Journal of Real Estate Finance and Economics, 14(3), 333-340.↩︎